The excitation operator approach to non-Markovian dynamics of quantum impurity 

models in the Kondo regime 
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We present a numerical method for studying the real time dynamics of a small interacting quantum 
system coupled to an infinite fermionic reservoir. By building an orthonormal basis in the operator 
space, we turn the Heisenberg equation of motion into a system of linear differential equations, 
which is then solved iteratively by constructing excitation operators. The application of our method 
depends on a layer structure in the operator space, which help us to turn an infinite linear system 
into a series of small systems. We apply the method to investigate the decoherence dynamics of 
quantum impurity models in the Kondo regime with a non-Markovian reservoir. Taking full account 
of environmental back-actions and electron-electron interactions, we find that the coexistence of 
the Kondo correlation and a non-Markovian reservoir induces coherence ringings, which will be 
suppressed by either driving the system away from the particle-hole symmetric point or changing 
the reservoir into a Markovian one. 

PACS numbers: 03.65.Yz, 02.60.Cb, 72.15.Qm, 73.23.-b 
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I. INTRODUCTION 

The decoherence of a small quantum system coupled 
to a fermionic bath has recently attracted much atten- 
tion [l|-[T7j. due to the fact that the fermionic bath 
manifests as an important source of decoherence in a 
wide range of electronic devices designed for solid-state 
quantum computers. In spite of considerable effort, a 
thorough understanding of the coherence dynamics of 
fermionic baths is still lack in the non-Markovian regime. 
The non-Markovian dynamics is difficult to address the- 
oretically, because the traditional Born-Markov approx- 
imation is invalid when the relaxation time of the en- 
vironment is comparably long and then the back-action 
of the environment plays an important role in the dy- 
namics of the system. To fully take into account the 
back-actions, the system-environment coupling must be 
treated in a non-perturbative way. In recent years, sev- 
eral non-perturbative approaches have been suggested 
to derive the master equation in the existence of strong 
back-actions [Iol - [l7j |. when the electron-electron interac- 
tion is absent or irrelevant to the non-Markovian dynam- 
ics. 

The non-interacting models, however, fail to incor- 
porate the physics in solid-state structures where the 
Coulomb interaction between electrons is greater than 
the electron kinetic energy. A well known paradigm is the 
Kondo effect, displayed in quantum dots in the Coulomb 
blockade regime. In the Kondo effect, the e-e interaction 
induces a strong correlation of electrons, which can only 
be understood from a many-particle point of view. Then 
it is obliged to study the interplay of correlation physics 
and non-Markovian dynamics. 

In this paper, we study the coherence dynamics of 



quantum dots in the Kondo regime coupled to a non- 
Markovian fermionic reservoir. The model is described 
by the Anderson impurity Hamiltonian, which can be 
written as 



H — Hs + Hb + H\ 



(1) 



is the sys- 



Here H s = €dJ2a%a d o 

tern Hamiltonian, where denotes the gate potential 
and U the Coulomb repulsive interaction. And Hb = 
Sfco- e kc\. a Cka is the Hamiltonian for a non-Markovian 
fermionic reservoir, which is set with a finite bandwidth 
and a sharp edge. The coupling Hamiltonian is given by 

Hv = J2ka V k (cQaCko + h.c), 
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To solve this problem, we develop the numerical exci- 
tation operator method on the basis of previous works by 
the author [TH, [l9[ . This method is designed for study- 
ing the real time dynamics of a strongly-correlated sys- 
tem driven out of equilibrium. As for quantum impurity 
models, it is distinguished from various approaches [20l — 
HH devised recently by the fact that both the Coulomb in- 
teraction and the system-environment coupling are dealt 
with in its full extent and at the same time the reservoir 
is set to be infinite. These are difficult to be fully realized 
in present approaches. 

The plan of the paper is the following. In Sec. II 
we introduce the excitation operator method. Its ap- 
plication in quantum impurity models is demonstrated 
in Sec. III. The results are discussed in Sec. IV. Espe- 
cially, we will discuss the intrinsic correlation between 
the non-Markovian dynamics and the Kondo physics. In 
Sec. V, we discuss the suppression of the non-Markovian 
dynamics. We conclude with a summary and discussion 
of our method and results in Sec. VI. 
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II. THE EXCITATION OPERATOR METHOD 

The excitation operator method is designed for solving 
the Heisenberg equation of motion: 



dO{t) 
dt 



i[H,6(t)], 



(2) 



where O is the observable that we are interested in. 

We choose an orthonormal basis {Oi} in the operator 
space which contains all the linear operators mapping the 
Hilbert space into itself. Any two basis operators satisfy 



A A) =5 t 



(3) 



where the bracket denotes the inner product between two 
operators and is generally defined as 



(diA 



i 



(4) 



Here A/" is the normalization factor. An arbitrary ob- 
servable can be decomposed into the linear combination 
of the basis operators. Then our target is to solve the 
Heisenberg equations of the basis operators. 

The Heisenberg equation is solved by constructing the 
excitation operators Ai satisfying the eigen equations: 



[H,A i ]=X i A i , 



(5) 



where H is the Hamiltonian of the system, and Aj the 
excitation energy of Ai. We suppose that the excitation 
operators are expressed by the basis operators as 



Ai — AijOj. 

3 



(6) 



The coefficients matrix A needs to be determined. We 
then calculate the commutators between the Hamiltonian 
and the basis operators 



[H,Oi] 



3 



(7) 



and obtain a matrix H. By substituting Eq. [S] and [7] into 
Eq. we find that A is in fact the unitary transformation 
to diagonalize the matrix T-L. 

In principle, the elements of T-L can be written as 



Hi, s = (di,[HA])- 



(8) 



By using the definition of the inner product and the fact 
that the Hamiltonian is self-adjoint, we prove that T-L 
must be a Hermitian matrix and diagonalizable. 

The solution of the Heisenberg equation can be ex- 
pressed by the coefficients A as 



Oi(t) 

Here we use the fact that A is unitary. 



(9) 



In practice, the dimension of W grows exponentially 
with the system size, so that directly diagonalizing it is 
impossible. However, there is a layer structure in the op- 
erator space, generated by the superoperator [H, ■], This 
indicates that we could change the problem of diagonal- 
izing T-L into the problem of diagonalizing a series of small 
matrices. 

We consider the evolution of the basis operator Oi in 
a small time interval r. The solution of the Heisenberg 
equation Oi(r) is mostly limited in a subspace of the 
whole operator space, generated by Oi and [H,Oi\. As 
t — } 0, we can calculate Oi(r) in this subspace, the di- 
mension of which is small. In other words, we express 
[H, Oi] as 



(10) 



where Hji is non-zero as j ^ i. Obviously, H is a sub- 
matrix of %. As t — > 0, the solution can be written as 



Um Oi(r) 

T— >0 



LAV 

3,i' 



f A 



? T A - iO 



(11) 



where Ai and A are the eigenvalues and the unitary ma- 
trix of T-L respectively. 

To calculate Oi(t) at a finite time, we divide the time 
t into N small intervals of length r = t/N, and have 



Oi(t) 



iHr 



-iHt 



(12) 



In each time interval, the evolution of the basis opera- 
tors is calculated according to Eq. [IT] As r — >• 0, the 
result will go to the solution of the Heisenberg equation. 
The point of this method is to utilize the layer struc- 
ture in the operator space. That is, the whole basis can 
be gradually generated from a single basis operator Oi 
by iteratively calculating the commutators between the 
Hamiltonian and the basis operators. But there is no 
such a structure in the Hilbert space. This is why we 
choose to solve the Heisenberg equation, instead of the 
Schrodinger equation. 

The number of basis operators which need to be stored 
in calculating Oi(Nr) increases exponentially as N in- 
creasing. This is a problem for numerical calculations. 
It can be solved by a truncation scheme: after obtain- 
ing Oi(Nr) one keeps only the M basis operators with 
the largest amplitudes. In this way, the number of the 
stored basis operators is fixed to be M and the compu- 
tation time increases linearly with N. Suitable values for 
the parameter M depends upon the model. It should be 
decided numerically by varying M. 
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III. THE MODEL AND THE ORTHONORMAL 
BASIS OF THE OPERATOR SPACE 

Our system consists of an impurity site, which is cou- 
pled, via particle-particle exchanges, to an electron reser- 
voir. The system plus environment is described by the 
Hamiltonian [Tj. To facilitate applying the excitation 
operator method, we re-express the Hamiltonian of the 
reservoir in real space in terms of an infinite chain: 



oo 



E 

(T.i— 1 



H B = ~9 V (cl a c i+ha + h.c), 



(13) 



where Ci a is the electron annihilation operator at site i. 
We take the size of the reservoir to be infinite. This 
avoids any coherence oscillation due to the finite-size ef- 
fects. The coupling Hamiltonian now becomes 



(14) 



where the system is coupled only to the first site of the 
reservoir. The system Hamiltonian keeps invariant as 



(15) 



We find an orthonormal basis in the operator space of 

3 

this model by transforming it into a spin-— chain by the 

Jordan- Wigner transformation. 

The model contains a series of sites. The dimension 
of the local Hilbert space at each site is four with the 
basis vectors 0), | t), I I) an d | t4-)- Keeping in mind 
that the basis operators are orthogonal to each other, we 
choose next sixteen 4x4 matrices as the local basis op- 
lOWl Wo° 0Wo° 
1 )' \ -1 J' ^ o a )' \ -o ( 

0lW0-lW0o Q \ H f0-o Q 
)> \o a ) am \o a 
the 1 denotes the two-dimensional identity matrix, and 
o a with a — x, y, z the three generators of SU(2) algebra, 
which are 



erators: 



Here 



1 

1 



-1 

1 



1 

-1 



(16) 



where L denotes the total number of sites and is taken as 
L — > oo. It is easy to prove that the operators in Eq. [17] 
satisfy the orthogonal relations. 

The Hamiltonian is real, so that % is a real symmetric 
matrix. And the diagonal elements of % are all zero. 
Because the basis operators are either symmetric or anti- 
symmetric, then the projection of [H,Oi] on Oi is zero. 

The Hamiltonian can be expressed in the basis oper- 
ators by the Jordan- Wigner transformation, in which a 
phase factor is attached to each site to produce the anti- 
commutative field operators. They are 



n 



o z 



<t 11 ^ \ -O z 

j<i v / J 



1 / o x 



i / Q y o 



and 



"4 ^ n 



2\0 o 1 . 1 2 o y 



o z 



(19) 



~o z 

I f -° Z 
2\o z 



1 ( 



2 \ o z 



, (20) 



where i, j = 0,1,2, - • • denote the sites and 



o z 
-o z 



the phase factor at site j. 

The hopping term in the Hamiltonian can then be ex- 
pressed as 



lfo x 
2^0 -o x 

lfo y 

2 I o - a y 



® 



1/0 1 
2^10 

1/0-1 
2 I 1 



o x 





v 


o x 


o y 


' 





o y 


o z " 






L 



i+1 



-o z 
o z 



(21) 



And the system Hamiltonian becomes 



The sixteen matrices form a complete basis of the local 
operator space. 

A basis operator can be expressed as the tensor prod- 
uct of the local operators: 



6 = Yl ®<Tj 



(17) 



where oi denotes the local operator at site i. Now we 
explicitly define the inner product as 



1 



(OuOj) :=_Tr[OjO 



(18) 



- U 2e d + U ( 1 

Hs = £d + 7 — [ -1 



2e d + U ( o z 



U ( o z 



o z J Q 4 V ~o z /n 

(22) 

After solving the Heisenberg equation, we need to cal- 
culate the expectation value of the basis operators with 
respect to the initial state. This is done by transforming 
the basis operators into Majorana operators defined as 



7 



C i(J ± Ci a ■ 



(23) 
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The sixteen local operators at site i are in one-to- 
one correspondence with next Majorana operators: 1, 

t+ t- 

t- 4- 

li li > 



4+ 

it > 7, 



4+ 4- 

7; 7i 
4+ 4- 
7i it 



t+ t- t+ 4+ t+ 4 
t+ t- 4+ t+ T- 4- t+ 4 

1 1 li li > / * li li ' / 2 



t-„,4+ 



l l 
+~4- 



7 



7 

the same 



and 7.J + 7.J it^lt • They are not exactly 
since the product of an odd number of Ma- 
jorana operators, such as contains phase factors at 
the sites j < i. However, we can design an iterative 
algorithm to transform a basis operator into a product 
of Majorana operators. The algorithm begins from the 
largest site where the local operator is not the identity, 
and sweeps the chain in the descending order. 

After the transformation, the expectation value is cal- 
culated by using the Wick's theorem. The contraction 
of a pair of Majorana operators at zero temperature is 
found to be 



w + ij + ) = -drir)= 6 ^> 



and 



(tTtT) = 



-2sin(|z - j|tt/2) 



\i-j\tt 



(24) 



(25) 



as \i — j\ is an odd number. 



IV. INTERACTION-INDUCED COHERENCE 
RINGING IN A NON-MARKOVIAN 
ENVIRONMENT 

We study the coherence dynamics of the system after 
its coupling to the reservoir is switched on at the time 
t = 0. The reduced density matrix of the system is ob- 
tained by calculating the expectation values of the sixteen 
local operators, and is formally written as J2i j Pn " I*} 0" I > 
where i,j = 1,2,3,4 and the corresponding states are 
|0), | t)j I I) an d | tl) respectively. Our method is distin- 
guished from the master equation approach by the fact 
that no approximation is made on solving the Heisenberg 
equation and the environmental back-actions are fully 
taken into account. 

We set the reservoir at zero temperature, avoiding the 
thermal fluctuation which would suppress the Kondo res- 
onance. The Fermi energy of the reservoir is set to be the 
energy zero. We employ the level-broadening at the im- 
purity site T, generally defined as T — V 2 /g [2^, as the 
energy unit. This is usually in studying the Anderson 
impurity model. And the time unit is set to be 1/r (the 
convention h = 1 is used throughout the paper) . 

At the particle-hole symmetric point, i.e., e ( i = — U/2, 
a large U provides a limit to the electron number of the 
system. The system is in the Kondo regime and can be 
described by a single spin. We suppose that its initial 
state is prepared as a superposition of the spin up and 
down states, i.e., a\ f) \). In the decoherence theory 
considering a Markovian environment, the coherence of 
the initial state will be lost in an exponential way after 
coupled to the reservoir. However, the real environment 




FIG. 1: Time evolution of the elements of the reduced density 
matrix at the symmetric point as U — 0. The results by the 
excitation operator method, represented by the black circles, 
are compared with the exact solution, represented by the vari- 
ous types of lines. The initial state is set to be (| T) + | 4-))- 
The coupling in the reservoir is set to be g — t. 



in the experiments is usually not Markovian, and the 
back-actions from the environment to the spin cannot be 
neglected. Here we consider a non-Markovian reservoir 
by setting the bandwidth of the reservoir to be compara- 
ble with the level-broadening at the impurity, i.e., g ~ Y . 

We first set the interaction U to zero to compare our re- 
sult with the exact solution, obtained by exact diagonal- 
ization of the single-particle eigenmodes. The elements of 
the reduced density matrix are shown in Fig. 1. The re- 
sult by the excitation operator method fits well with the 
exact solution, until the density matrix has relaxed to its 
equilibrium value. This proves that our method is a pow- 
erful tool in studying the real time dynamics of a quan- 
tum system coupled to a non-interacting reservoir. The 
errors can be controlled by letting r — > and M — > oo. 
It provides a reliable way of understanding the dynamics 
of decoherence, especially in strongly-correlated systems, 
where no analytical method is available. 

As the e-e interaction is absent, the off-diagonal ele- 
ment, i.e., the coefficient of the term | t)(i \> decays ex- 
ponentially, as predicted by the decoherence theory. This 
is the feature of a Markovian dynamics. The back-actions 
of the environment are sufficiently suppressed. However, 
it is not the case as U ^ T (see the top panel in Fig. 2). 
In the existence of a strong interaction, the exponential 
decays are replaced by oscillations. And the intermedi- 
ate quasi-steady regimes are observed. The decoherence 
time significantly increases as U increasing. We then an- 
alyze the time evolution of the von Neumann entropy at 
different U (see the bottom panel in Fig. 2). As is well 
known, the equilibrium value of the entropy is 2 In 2 [39j . 
As U — 0, the entropy increases monotonically from zero 
towards its equilibrium value, corresponding to the expo- 
nential decay of the off-diagonal elements. But as U » T, 
we find strong oscillations in the entropy, which is the 
signal of non-Markovian dynamics. The coherence in the 
initial state is lost and recovered repeatedly, similar to 
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Time(l/r) 

FIG. 2: The off-diagonal element and the von Neumann en- 
tropy of the density matrix as a function of time at different 
U (top: off-diagonal element, bottom: entropy). We choose 
the initial state of the system to be t) + 4-)- The 

entropy for a different state y\\ t) + \J~\ \ i) is a l so studied 
at U — 16r, which is represented by the black circles for a 
comparison. 

the spin echo effect. However, in our model, the purifi- 
cation of states arises naturally from the e-e interaction 
and no external driving field is needed as in the spin echo 
or dynamical decoupling technique (40| . This provides a 
new perspective in protecting the quantum state. 

The coherence ringing is an effect induced by the e-e 
interactions. It must be distinguished from the oscil- 
lations of coherence observed in the non-Markovian en- 
vironments before [ill . [l2j , where the e-e interaction is 
absent. Without interactions, the electrons move inde- 
pendently, and the non-Markovian dynamics can be un- 
derstood in the single particle picture. However, in the 
existence of a strong interaction, the single-particle pic- 
ture breaks down due to the correlations between elec- 
trons. As U ~3> T, the dissipation process is controlled by 
the Kondo correlation. The steady state as the time goes 
to infinity is a spin singlet. The correlation between the 
spin in the system and the spins in the reservoir is built 
in course of time, accompanied by the loss of coherence 
in the system. The non-Markovian coherence dynamics 
is in fact related to the dynamics of spin correlations in 
a Kondo model. 

The initial state is found to be indifferent to the dy- 
namics of decoherence. As an example, we choose two 
different initial states, which are A=\ t) + ^1 4-) an d 

f) + 4-), and plot the corresponding entropies 




1.4 1.6 1.8 2 2.2 2.4 2.6 2.8 
ln(U/T) 



FIG. 3: In this figure, we plot the period of coherence ringing 
obtained from the numerical result, represented by black cir- 
cles, as a function of U. The data is fitted to a linear function, 
represented by the solid line. 

as a function of time in Fig. 2. The results are exactly the 
same. This reflects the spin-flip symmetry in the Hamil- 
tonian. The correlation between the non-Markovian dy- 
namics and the Kondo physics is universal for the initial 
state a\ 7) + P\ I)- 

The period of the coherence ringing is obtained from 
the numerical result, and plotted as a function of U in 
Fig. 3. In a large regime of U, the oscillation period is 
found to be proportional to In (U/T). As the interaction 
strength increasing, the period decreases in a logarithmic 
way. That the oscillation period depends on U can be 
understood by studying the energy levels of the system. 
We emphasize that the coherence ringing happens as the 
system is at the particle-hole symmetric point, i.e., = 
—U/2. At the symmetric point, the Hamiltonian of the 
system changes into 

/ 1 \ 

( S -o 1 A S ■ « 

V 1 / 

We see that the energy levels of the system are degener- 
ate. The ground level is two-fold degenerate, containing 
the spin up and down states. It is separated by a gap 
of U/2 from the excited level, which is also two-fold de- 
generate, containing the vaccum state and the doubly- 
occupied state. An energy gap of U /2 protects the sub- 
Hilbert space containing the states | 7) and | 4-) , and then 
protects the quantum coherence in the initial state. The 
gap is critical to the appearance of coherence ringing, 
which is obvious only as the gap is large. 

V. SUPPRESSION OF THE COHERENCE 
RINGING 

We attribute the coherence ringing to the coexistence 
of the Kondo correlation and the non-Markovian reser- 
voir. Then it should disappear if any of the two condi- 
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FIG. 4: The time evolution of the entropy at different e d . 
The interaction is set to be U = 8r. Then the particle-hole 
symmetric point is at e d — —4F. 
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FIG. 5: The time evolution of the entropy at the particle-hole 
symmetric point at different g. The interaction is set to be 
U — 8r. The bandwidth of the reservoir g — > oo corresponds 
to the Markovian limit. 



tions is broken. This is verified by the numerical results 
(see Fig. 4 and 5). 

In Fig. 4, we plot the time evolution of the entropy at 
different gate potentials. As the system is away from the 
particle-hole symmetric point, the coherence ringing is 
suppressed. This is due to the suppression of the Kondo 
resonance as the system is depleted or doubly-occupied. 
The suppression of coherence ringing can also be under- 
stood by the splitting of the excited level. The coherence 
ringing is distinguished from a simple Rabi oscillation be- 
cause there are totally four levels in the system. By driv- 
ing the system away from the symmetric point, we break 
the degeneracy at the excited level, which splits into the 
vacuum level and the doubly-occupied level. This is re- 
lated to the disappear of the coherence ringing. Accord- 
ing to our knowledge, it is the first time to find that the 
coherence in the ground state depends on the degeneracy 
at the excited level. 

The e-e interaction induces a coherence ringing only 
if the environment is non-Markovian, i.e., g ~ T. In 



Fig. 5, we show the entropy functions at different g, the 
bandwidth of the reservoir. In the case of g = T, the non- 
Markovian dynamics is significant. The coupling between 
the neighbor sites in the reservoir is as same as the cou- 
pling between the system and the reservoir. The back- 
action is strong since the relaxation time in the system 
is comparable with that in the reservoir. As the inter- 
action is much larger than the bandwidth, i.e., U 3> g, 
the coherence ringing appears. If we keep the interaction 
invariant, at the same time increasing the bandwidth <?, 
the coherence ringing is suppressed. As g ~ U 3> T, the 
oscillation is totally destroyed. In the limit of an infi- 
nite band, i.e., the Markovian limit, the entropy function 
recover the feature in the non-interacting model. That 
is, it increases monotonically towards the steady value: 
2 In 2. But the relaxation time is now controlled by the 
interaction U, instead of the impurity level- width T. We 
see the disappear of the coherence ringing as the reservoir 
changes gradually to the Markovian limit. 



VI. CONCLUSIONS 

We have presented the numerical excitation operator 
method to coherence dynamics of an interacting quan- 
tum system coupled to a fermionic reservoir. Compared 
with the present analytical approaches, it takes full ac- 
count of the Coulomb interaction between electrons and 
the system-environment coupling, and then provides new 
information on the interplay of electron-electron correla- 
tions and environmental back-actions. At the same time, 
our method takes into account an infinite reservoir by uti- 
lizing the layer structure in the operator space, and then 
avoids the finite-size effects which threaten the present 
numerical methods. We have applied this method to an 
interacting quantum dot coupled to a fermionic reservoir, 
discovering the coherence ringing induced by the e-e in- 
teraction in the Kondo regime. The coherence ringing 
is a many-body effect and can only be observed in the 
presence of both Kondo resonance and non-Markovian 
reservoir. It will be suppressed as the system is away 
from the particle-hole symmetric point or the reservoir 
changes into the Markovian limit. 

Although we concentrate in this paper the dynamics of 
decoherence in quantum impurity models. The method 
that we presented can be applied to investigate the real 
time dynamics in a wide range of models describing a 
quantum system coupled to spin, fermionic and bosonic 
reservoirs. 
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